Antiferromagnetism in doped anisotropic two-dimensional spin-Peierls systems 
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We study the formation of antiferromagnetic correlations induced by impurity doping in anisotropic 
two-dimensional spin-Peierls systems. Using a mean-field approximation to deal with the inter-chain 
magnetic coupling, the intra-chain correlations are treated exactly by numerical techniques. The 
magnetic coupling between impurities is computed for both adiabatic and dynamical lattices and is 
shown to have an alternating sign as a function of the impurity-impurity distance, hence suppressing 
magnetic frustration. An effective model based on our numerical results supports the coexistence of 
antiferromagnetism and dimerization in this system. 

PACS: 75.10 Jm, 75.40.Mg, 75.50.Ee, 64.70.Kb 



General interest for spin-Peierls (SP) systems was re- 
cently renewed by the discovery of CuGeOaEI, the first 
inorganic SP material. The SP transition is character- 
ized by a freezing of the spin fluctuations below an en- 
ergy scale given by the spin gap AjS accompanied by 
a simultaneous lattice dimerizationB. Rich phase di- 
agrams have been obtained experimentally upop,.-dop- 
ing this compound with non-magnetic impuritiesau. In 
site-substituted systems such as (Cui_ ;c M a: )Ge03, where 
M=Zn (Ref. ||) or Mg (Ref. |), long range antiferro- 
magnetic (AF) order is stabilized at low temperature 
while the dimerization still persists (D-AF phase). In 
Mg-doped compounds, for impurity concentrations larger 
than a critical value (x c — 0.02), a first order transition 
occurs between the D-AF phase and a uniform AF (U- 
AF) phase where the dimerization disappearsa The co- 
existence of the two types of order in the D-AF phase is 
an intriguing phenomenon since lattice dimerization fa- 
vors the formation of spin singlets on the bonds while low 
energy spin fluctuations exist in an AF phase. 

Theoretically, the effect of impurity doging in SP 
terns was considered for, fixed-dimerizecH, adiabaticE 
and quantum-dynamicala lattices. A single nonmagnetic 
impurity releases a soliton in the chain which can be 
viewed as a kink in the lattice distortion. In the absence 
of interchain couplings, such an excitation can freely 
propagate away from the impurity. On the other hand, 
the interchain elastic coupling K±_ was shown to produce 
confinement within some distance from the impurityou. 

For a finite impurity concentration, the coexistence be- 
tween SP and AF orders has been previously discussed 
either considering randomly distributed domain walls in 
a XX chaina or assuming-small fluctuations of the mag- 
netic exchange constantaij. Despite their success to de- 
scribe some experimental results, these models are rather 
limited since they do not take into account the micro- 
scopic origin of the soliton formation nor the interchain 
couplings. In this paper, a realistic microscopic model 



with interchain magnetic and elastic couplings is consid- 
ered to describe the formation of a region with AF corre- 
lations in the vicinity of each impurity and which allows 
an estimation of the effective interaction between impu- 
rities in the two-dimensional (2D) system. Thus, we are 
able to construct and study an effective model in order 
to understand the effects of a finite impurity concentra- 
tion.B 

In a first step, the spin-phonon coupling is treated 
in the adiabatic approximation. The Hamiltonian TL = 

Hmag + Hel IS: 

Hmag = J\\ VJ(1 "b ^i,a) Sj, ' Si+i, a + J± S Si a ■ S^J , 



i,{a,b) 
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where a is a chain index and i labels the sites along 
the chains. Atomic displacements are only considered 
along the chain direction, i5^ a being here a classical vari- 
able related to the change of the bond length between 
sites (i,a) and (i + I, a). The magnetic part includes a 
magnetoelastic coupling Jii (hereafter set to unity) and 
an exchange interaction J± connecting nearest neigh- 
bor (NN) chains. We eventually include in our model 
a next NN exchange interaction along the chajp-, whose 
relevance for CuGeOa has been emphasizedEJH. TC e i 
is the elastic energy. The interchain elastic interaction 
(K±) is limited to NN chains. Stability of the lattice 
implies K» > 2|2fj_|. Typical values of the parameters 
for CuGe0 3 arc Jj_ ~ 0.1 (Ref. ||) and K x /K\\ - 0.2 
(Ref. If). 

In order to study numerically model (||), we treat ex- 
actly the single chain problem using exact diagonaliza- 
tion (ED) or Quantum Monte Carlo (QMC) methods, 
while the interchain magnetic coupling is treated in a 
self-consistent mean-field (MF) approximation. This is 
an standard procedure to include interchain couplings in 
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the study of quasi-one-.dimensional systemslll Moreover, 
Inagaki and FukuyamatZI have used a similar MF approx- 
imation to treat the interchain coupling in the bosonized 
version of (^) within a self-consistent harmonic approxi- 
mation. Thus, in our procedure, the interchain magnetic 
coupling is replaced by its MF form: 



H 



MF 



i, (a, 6) 



SZaW-iSZaKSZb)}- (2) 



By extending a similar approach previously applied to 
one-dimensional (ID) chainst§B to the case of the 2D 
lattice, a sweep is performed in the transverse direction, 
i.e. a — ► a + 1. For each chain a, we compute the MF 
values (S* a ) and the classical variables {S ita } by energy 
minimization, which is achieved by solving iteratively the 
equations 

5%a = —{J\\{Si,a ■ S i+l,a) + Kx(Si, a +i + 5i, a -x)} / K\\. (3) 

Then, these new values of the AF and SP order parame- 
ters enter as input for the chain a + 1. This procedure is 
iterated until convergence is reached. In this way, we can 
study numerically finite clusters consisting of N coupled 
chains with L sites where, typically, N x L = 12 x 18 in 
ED and iV X L = 6 X 40 in QMC, with toroidal boundary 
conditions. 

A similar MF approach can be adapted to study a 
model equivalent to (|]) but with quantum phonon de- 
grees of freedom. In this case, phonon operators b\ a and 
bi :a are introduced on each bond and the displacements 

Si y a become g{b\ a + &i, a ), where g is the magnetoelastic 
constant. Then, the classical elastic term Ti. e i is replaced 
by its quantum version, 

n ph = nY / {bl^, a + r(b tM + blj(b lM+1 + bl a+1 )} , (4) 



where T = Kx/(2Kn), and the phonon frequency f2 is 
related to g by ft — 2g 2 K\\ . The adiabatic limit ([j]) is 
recovered when f2 — > (requiring g — + also) . Similarly 
to the interchain magnetic term, the interchain elastic 
term of is then treated in mean-field by introducing 
a lattice order parameter 8fa = 9$i a + h,a}- Then, the 
term (0) is replaced by: 



r 



(5) 



Note that in this case it is not necessary to solve an equa- 
tion similar to (||) . To diagonalize the single chain Hamil- 
tonian with L < 8, a Lanczos algorithm is used. The 
phononic degrees of freedom are treated-,;&4thin a varia- 
tional formalism previously introduced 
inelastic neutron scattering experiment, 




Note that 
on CuGeC>3 



reveal a rather large phonon frequency Q/Jii ~ 2 sug- 
gesting large lattice quantum effects in this material. 



As a preliminary study, we apply the MF procedure to 
the case of an homogeneous system without impurities. 
In models ([[]) and (|4|), J± is expected to stabilize the AF 
state while a small (or large magnetoelastic coupling 
g) tends to favour SP order. For each value of the cou- 
plings J± and K (where K — K» — 2\K±\ is the relevant 
parameter in the SP phase) we obtain the ground state 
without imposing any restriction on the MF parameters. 
We found only two different phases, the SP phase where 
(S z ) = and Si.a ^ and the antiferromagnetic (AF) 
phase with (S z ) ^ and 8i :a — 0. Then, the phase di- 
agram in the K — J± plane, can be obtained in a more 
efficient way by a direct comparison of the energies of 
the Neel antiferromagnetic phase, and of the uniforml 
dimerized phase. The phase diagram shown in Fig 
exhibits a transition line between AF and SP phases. 
In the adiabatic case, this line could be fitted by a law 
Jx = £+B with A = 0.3656 and B = -0.06 (this arti- 
ficial small negative value may be a consequence of small 
finite size effects, see e.g., Ref. |l^). A phase boundary 
with a foriiL J± = £ has been predicted by Inagaki and 
Fukuyamali3. However their bosonized approach do not 
fix unambiguously the value of A. 
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FIG. 1. Phase diagram of the pure system as a function of 
K — K\\ — 2\Kx\ and J±. The full lines correspond to the 
adiabatic results obtained by ED for chains with up to L = 18 
sites, q = (circles) and a = 0.36 (squares). The dashed 
(dotted) line corresponds to the calculation with quantum 
phonons for L = 8, K± = 0.2 a = 0, and fi = 1 (fi = 2) 

In the case of the adiabatic calculation finite size effects 
were shown to be small for K < 2. On the other hand, for 
the dynamical lattice, the calculation is reliable only for 
larger lattice couplings (i.e. smaller values of K) due to 
stronger finite size effects. As expected, for very small K, 
lattice quantum fluctuations are less effective in dimeriz- 
ing the chain than the adiabatic lattice. Then, the phase 
boundary obtained with quantum phonons is located be- 
low the adiabatic one. This tendency becomes clear as 51 
is increased, as shown in the figure. On the other hand, 
it has been suggested that dynamical phonons induce an 
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effective magnetic frustrationES. This frustration which 
becomes relatively important for larger K destabilizes 
the AF phase thus moving the phase boundary upwards 
as seen in the figure. Consistently with this behaviour, 
if a next NN exchange term is included in the Hamil- 
tonian, in the adiabatic approximation, the SP phase is 
more stable and the phase boundary is located above the 
corresponding curve for a — 0. This behavior is shown 
in Fig. |l| for the realistic value of a = 0.36 obtained 
for CuGeC>3 (Ref. |l^), where a is the value of the next 
NN exchange coupling constant in units of J|| . We have 
checked that the set of realistic parameters K ss 20 and 
K±, J±_ mentioned above corresponds to a point in the 
SP phase. 

To start our analysis of impurity doping, we consider 
a single impurity in order to investigate the appearance 
of AF correlations in the SP phase. As mentioned above, 
the impurity releases in n the chain a topological spin- 
1/2 solitonic excitatiorou characterized by a change of 
parity of the dimerization order which occurs in a fi- 
nite region of longitudinal size £m given by the soliton 
width. The local magnetization on each chain a can Jbe 
decomposed into uniform and staggered componentsE3, 
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In fact, the excess uni- 



form component S* ol = ±| and the soliton, characterized 
by a broad maximum of M?^ 8 , remain confined in the 
chain with the impurity. However, as seen in Fig. ||(a), 
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FIG. 2. Typical patterns of the local staggered magne- 
tization M^ s along the chains. The transverse coordinate 
a is indicated for each curve (only a few chains of the 2D 
cluster are shown). Closed circles correspond to the chains 
with an impurity, (a-b) QMC results for a 8 x 40 cluster with 
Ku = 1.9, Kj_ = 0.2 and J± — 0.1; (a) single impurity case 
(io = 21); (b) 2 impurities on the same chain at a distance 
Ai = L/2 = 20 (ii = 1 and i 2 = 21). (c-d) ED results for two 
impurities located on next NN chains at a distance of Ai = 9 
along the chain (as indicated by the arrows) on a 12 x 18 
cluster for K\\ = 2.4, K± = 0.2 and J± = 0.1. Singlet (d) and 
triplet (c) configurations are shown. 



the interchain magnetic coupling J±_ generates a large 
staggered component with the same parity, i.e. Af l st ag 
keeping the same sign, in the neighboring chains. Si- 
multaneously, the amplitude of the SP dimerization is 
significantly suppressed compared to the bulk value i.e. 
far away from the impurity. Large AF correlations can 
be seen up to more than four chains away from the impu- 
rity chain for magnetic couplings as small as Ji = 0.1, in 
particular in the vicinity of the SP — > AF transition line 
of Fig. [p. The transverse range of the AF 'polarization 
cloud' around the impurity increases strongly with the 
transverse coupling J±. 

A crucial feature of the polarization surrounding the 
impurity-soliton area is that the sign of is unam- 

biguously fixed by the orientation (5* ol = ±5) of the soli- 
ton and by the position («Oj a o) of the impurity in such 
a way that sign^M^} = sign{S^ ol } (-l) ao+io+1 . This 
fact can be simply understood in the strong dimerization 
limit (6 — > 1) where the introduction of the impurity on 
a given site releases a spin- 1/2 on one of its neighboring 
sites by breaking a singlet bond. For smaller lattice cou- 
pling, the excess spin can effectively hop from site to site 
on the same sublattice (due to the underlying dimeriza- 
tion), hence producing AF correlations with the parity 
defined above. 

Let us now consider two impurities introduced simul- 
taneously on two sites («i,Oi) and (121^2) of different 
chains (ai 7^ a%). When the two polarization clouds as- 
sociated to each soliton-impurity 'pair' start to overlap, 
one expects their interaction to depend on the relative 
orientation of the two solitons. As seen in Figs. ||(c,d), 
quite different patterns correspond to the singlet and 
triplet arrangements of the two spin- 1/2 solitons. As 
confirmed by our calculations, the lowest energy is always 
obtained for a spin state which leads to the same parity 
of the staggered magnetization associated to each impu- 
rity, i.e. which avoids completely magnetic frustration. 
The simple argument developed above for a single impu- 
rity then suggests that a triplet S = 1 (singlet 5 = 0) 
configuration is favored when the two impurities are lo- 
cated on the same sublattice (opposite sublattices). It 
is then appropriate to define an effective magnetic cou- 
pling between the AF clouds associated to each impurity 
by J e s = Eg=i — Es=o- For a wide range of parameters 
leading to a SP state in the bulk, we have numerically 
found that J c g is ferromagnetic (F) if the two impuri- 
ties belong to the same sublattice and antifcrromagnetic 
in the opposite case. This implies that the coupling be- 
tween the two local moments associated to the impurities 
is either F or AF in such a way that no frustration oc- 
curs. The magnetic coupling, for physical values of the 
parameters, can be fairly extended in space as seen in 
Fig. |](a). Its range is directly controlled by the overlap 
of the polarization clouds. It follows roughly a behavior 
like 

Jeff w J (-l) Aa+Al+1 exp (-C x p) exp {-C\A, (6) 

U 4|| 
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where Jo ~ J_l, C a are of order unity and £ Q ~ J a / A$. 

To get an insight on how finite size effects might affect 
our results we have increased our cluster size (12 x 18) in 
both directions. The change of the transversal dimension 
has very little effect because the polarization clouds are 
almost independent when the impurities are separated 
by more than four chains (see the very small values of 
J c ff for Aa = 4 in Fig. ||(a)). On the other hand, by 
increasing the length of each chain, some changes occur 
but only when the impurities are located at the largest 
distances. However, the exponential decay of the effective 
interactions with distance (see below) is not qualitative 
changed as it can be seen in |^(a)). Therefore, we expect 
that the numerical values of the fitting parameters would 
not be much affected by finite size effects leaving the 
overall behaviour essentially unchanged. In particular, 
we believe that the presence of long range AF order in the 
effective model (see discussion below) is a robust feature 
not affected by finite size effects. 
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in agreement with previous workfcl and the triplet excita- 
tion energy remains large (~ As). Then, one can expect 
that for larger chains, when the formation of soliton be- 
comes favourable, the effective interaction between them 
will be AF and their magnitude will be of the order of 
Ag. In summary, one can assume that the effective in- 
teraction between impurities on the same chain has the 
same form of Eq. (^|) with Aa = and Jo being now 
w As- This form is similar to the one adopted in Ref. ^ 
for impurities in a single chain except that these authors 
do not consider the sublattice sign alternation. Neverthe- 
less, it should be noticed that this AF interaction should 
also decay for very large Ai. 

When lattice quantum fluctuations are introduced 
(f2 > 0), the qualitative properties of the effective in- 
teraction J e ff are preserved. Consistently with the rela- 
tively larger stability of the AF phase in the small K 
region shown in Fig. [I] with respect to the adiabatic 
case, we have found that lattice dynamics lead to an in- 
crease of the size of the AF cloud associated to each soli- 
ton. Therefore, the magnitude of the magnetic coupling 
J G ff increases with the phonon frequency SI as shown in 
Fig. |(b). 

The final part of our study is the analysis of a sim- 
ple effective two-dimensional spin- 1/2 Heisenberg model 
between impurities with a long range interaction given 
by Eq. (^). The 'bare' parameters are the same as in 
Fig. || and the parameters of the expression (Q) have been 
obtained by fitting the curves shown in that figure and 
similar data for the case of Aa = 0. In the direction per- 
pendicular to the chains we have neglected the effective 
interactions beyond a distance Aa = 5. We have also 
assumed that even segments are associated to a soliton- 
antisoliton pair a 
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FIG. 3. Absolute value of the effective magnetic coupling 
between two impurities located on different chains separated 
by Aa; (a) vs. the separation Ai = «2 — ii along the chains 
for an adiabatic lattice (SI — 0) and for Kit = 2.4, K± — 0.2 
and J± = 0.1; (b) vs. the frequency SI for Ku = 1, K± — 0.2 
and Jj_ = 0.1. The sizes TV x L of the clusters used in the 
calculations are indicated on the plots. 

When two impurities are introduced on the same chain 
(ai = 02) two cases have to be distinguished. If the impu- 
rities are located on the same sublattice a similar behav- 
ior is observed as described above (compare Fig. ||(b) and 
Fig. ^(d)), i.e the effective interaction is /erromagnetic. 
However, the magnitude of |J e //| is ~ 0.4 (very slowly 
decaying as Ai increases) for the parameters of Fig. ||(a), 
i.e. much larger than the values corresponding to impu- 
rities in different chains, ff the impurities belong to dif- 
ferent sublattices then a chain with even number of sites 
is cut into two segments with even number of sites each. 
In the lowest energy configuration (S = 0) no soliton- 
antisoliton pair was observed for separations Ai up to 20 
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FIG. 4. Square root of the staggered magnetization (de- 
fined in the text) as a function of the impurity doping for 
several 2D clusters. Full circles correspond to, from top to 
bottom, 12 x 12, 16 x 16 ,■ • • , 40 x 40 clusters. The open 
circles are the extrapolations to the bulk limit. 



4 



A given number of spin- 1/2 impurities 4 < iVj mp < 
16 is thrown at random on systems of coupled chains of 
sizes up to 40 x 40. Then, the staggered magnetization 
M stag = (l/JV?)(E <) „(-l) <+0 'S? ) a J wh ere N s = N x L, 
is computed and averaged over, typically, 12,000-16,000 
random realizations. The square root of this quantity is 
shown in Fig. |[ By extrapolating to the bulk limit for a 
fixed impurity doping using a polynomial in (l/y/N a ) we 
found that Af sta g is finite, implying long range AF order, 
and slowly decreasing as x goes to zero. This behaviour 
is consistent with experimental results! suggesting that 
-Mstag decays exponentially to zero as x — > 0. 

In conclusion, spin-1/2 solitons released in 2D 
anisotropic SP systems by the introduction of impurities 
were shown to experience spatially extended F or AF ex- 
change interactions depending on their relative positions. 
These exchange interactions coexisting with the SP or- 
der are calculated from realistic microscopic models and 
used to construct a simple effective model which in turn 
enables us to show the establishment of long range AF 
order and to compute the AF order parameter as a func- 
tion of the impurity doping. 

This work was supported in part by the ECOS-SECyT 
A97E05 program. We thank IDRIS (Orsay, France) and 
Florida State University for using their supercomputer 
facilities. 
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